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Abstract 

The Christensen-Egan algorithm is extended and generalized to efficiently eval- 
uate new spin foam vertex amplitudes proposed by Engle, Pereira & Rovelli 
and Freidel & Krasnov, with or without (factored) boundary states. A concrete 
pragmatic proposal is made for comparing the different models using uniform 
methodologies, applicable to the behavior of large spin asymptotics and of ex- 
pectation values of specific semiclassical observables. The asymptotics of the new 
models exhibit non-oscillatory, power-law decay similar to that of the Barrett- 
Crane model, though with different exponents. Also, an analysis of the semi- 
classical wave packet propagation problem indicates that the Magliaro, Rovelli 
and Perini's conjecture of good semiclassical behavior of the new models does 
not hold for generic factored states, which neglect spin-spin correlations. 

PACS numbers; 04.60.Pp 



1 Introduction 



Spin foam models are an attempt to produce a theory of quantum gravity start- 
ing from a discrete, path integral-like approach. They were first defined a decade 
ago (6j[l0]. More recently, we have seen significant progress toward extraction of their 
semiclassical behavior and its favorable comparison to the expected weak field limit of 
gravity, starting with Rovelh and collaborators' calculation of the graviton propaga- 
Unfortunately, further calculations have revealed that the standard spin 
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foam model due to Barrett and Crane produced incorrect results for some of the prop- 
agator matrix elements (3[|4]. This result has motivated several proposals to replace 
the Barrett-Crane (BC) spin foam vertex amplitude 10 for quantum gravity. The 



first proposal, by Engle, Pereira and Rovelli (EPR) |20 21 , aimed also to identify the 
spin foam boundary state space with that of loop quantum gravity spin networks; this 
model is also referred to as the "flipped" vertex model. Another proposal, by Livine 
and Speziale 27 28 , used S'/7(2)-coherent states to define the spin foam amplitudes 
and reproduced the EPR proposal up to an edge normalization factor. Finally, a paper 



by Freidel and Krasnov [18}|22] , suggested that the EPR model corresponds to a topo- 
logical theory related to gravity and proposed a generalization thereof corresponding 
to gravity itself (the FK model). The present paper, along with most previous work, 
concerns only the Riemannian signature models of gravity. 



Having been defined, the new models must be tested to see whether their semi- 
classical behavior is an improvement over the BC model. The standard tests involve 
examining the asymptotics of their vertex amplitudes and checking the semiclassical 
behavior of observables. So far, two test problems involving observables have been 
proposed: semiclassical wave packet propagation j29], and evaluation of the graviton 



2-point function 12 17 34 . Both problems require the computation of large sums, 
where the spin foam vertex amplitude is contracted with a suitably defined boundary 
state. These computations, while important for extracting the physical content of the 
new spin foam models, have so far not been tractable. 

This paper describes efficient numerical algorithms, based on the existing Chris- 
tensen-Egan (CE) algorithm for the BC model, to evaluate the new spin foam vertex 
amplitudes, both with fixed boundary spins and contracted with a boundary state. 
The accessible boundary states are restricted to the large (though with important 
limitations) class of so-called factored boundary states. Calculation of the new vertex 
amplitude asymptotics shows that they are dominated by the degenerate spin foam 
sector, as is the BC model. Also, the numerical simulation of semiclassical wave packet 
propagation shows that, for factored boundary states, the shape of the propagated 
wave packet does not agree with the desired semiclassical result. The state of the 



good semiclassical behavior conjecture 29 for states that do not neglect spin-spin 



correlations (unlike factored states) remains open. 

Section [2] briefly introduces spin foam models and presents the three models de- 
scribed above in a unified framework, also incorporating spin foams with boundary. 
Section [3] describes the proposed calculations of observables and asymptotics. The 
class of factored boundary states is tentatively defined in section 2.3 and in full detail 
in section [5] In section |4] explicit formulas are given for new spin networks arising 
in the evaluation of spin foam model vertex amplitudes. Section [5] presents efficient 
numerical algorithm for computation with the new models and section [6] shows the 
results of their application to some of the problems discussed in section|3] Finally, 
section [7] concludes with a discussion of the results and future work. 



2 BF theory and spin foams 

The new spin foam models of gravity may be presented in a way similar to the orig- 



inal BC model. Following Freidel and Krasnov 22 , we define them within a unified 



framework. See also the more recent paper 18 . 

The starting point is BF theory. It is a 4-dimensional field theory with two fields: 
a gauge connection 1-form A and an auxiliary 2-form B. The action is given by 

J tiiBAF], (1) 

where F — dA is the curvature of the connection. If the gauge group is taken to be 
Spin(4), the double cover of SO{A) and a constraint is imposed, ensuring simplicitjj^of 
the B 2-form, this theory becomes equivalent to the Plebanski formulation of general 
relativity in Riemannian signature [10| . 

BF theory is in a sense topological. Particularly, its underlying manifold may be 
freely changed from a smooth one to a discretized (piecewise linear) one, as long as the 



^Simplicity means that there exists a 1-form e such that _B = *(e A e), with * denoting the Hodge 
dual. 
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topology remains the same. Moreover, for BF theory, quantization and discretization 
commute [s]. Spin foam models aim to reproduce gravity by heuristically imposing 
simplicity constraints on BF theory after the discretization and quantization steps 
have been performed [T] . The connection between spin foams and gravity is motivated 
by results from loop quantum gravity [33] . 

Consider BF theory defined on a simplicial complex (possibly with boundary), 
also referred to as a triangulation. It is convenient to introduce the dual 2- complex. 
Each 4-simplex is identified with a dual vertex, each tetrahedron is identified with a 
dual edge, and each triangle is identified with a dual face. Discretizing the A and B 
fields and integrating out the B field, the theory's path integral yields the following 
expression for its partition function: 

- y dAdBe'^ " = j dg\{5{g}), (2) 

where the connection A has been replaced by group elements g associated to every 
dual edge, and gj represents the holonomy around a dual face. This form of the 
partition function manifestly shows that only flat connections (with trivial holonomies) 



contribute to the BF theory path integral. See 31 for details 



The (5-functions can be expanded in terms of gauge group characters and the group 



integrations can be performed at each dual edge 31 . What remains is a discrete sum 
of the form 

A{F) = \{Af{F)\{A,{F)\{A,{F), (3) 

f e V 

Z = Y.A{F), (4) 

F 

where F ranges over all spin foams (defined below), while /, e, and v range respectively 
over dual faces, dual edges, and dual vertices. In this context, a spin foam is a labelling 
of the dual faces of the triangulation by irreducible representations of the gauge group. 
These representation labels come from the character expansion described above. This 
definition of spin foams will have to be augmented with extra labels for the purpose 
of introducing the new models. Thus, discrete BF theory yields a spin foam model. 
In general, a spin foam model is defined by they way it labels spin foams and by the 
amplitudes it gives to them through (|3| and specific choices for Af{F), A(,{F) and 
Ay{F). Each of the amplitudes for dual cells may depend on its own label and the 
labels of adjacent cells. 

Irreducible representations of Spin(4) = SU{2) x SU{2) are labelled by a pair 
of integers j — where each j is a spiij^ corresponding to an irreducible 

representation of SU{2). Hence forth, all representation labels will be referred to as 
spins, unless otherwise specified. Bold face letters will specifically represent Spin(4) 
spins. 

For pure BF theory, face amplitudes are determined by the character expansion 
of (5-functions and are given by the dimension of the irrep j y labelling a given face 

Af{F) - dimjy - {jj + l)(j+ + 1). (5) 



Technically, a twice-spin, since it does not take on half-integral values. 



3 



Edge and vertex amplitudes are determined by evaluating the group integrals in equa- 
tion ([2]). The basic identities we use are 




(-1) 



(6) 



The above graphical notation requires some explanation. See 31 and the Appendix 
of 
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for full details. Briefly, a solid line represents a tensor index, while any object 
with lines attached to it is a tensor, with each line representing a tensor index. A 
vertical strand with a circle represents a matrix element of a particular representa- 
tion of the gauge group. Concatenation of lines corresponds to index contraction; 
particularly concatenating strands with circles corresponds to matrix multiplication. 
Juxtaposition of tensorial objects corresponds to their Kroncckcr (tensor) product. 
These graphical representations of SU{2) or Spin(4) tensor contractions will be re- 
ferred to as spin networks. Most of the representation labels and basis indices have 
been omitted for conciseness. Instead, some spins will be marked as collective labels 
with an asterisk. Their expanded meaning should be clear from context. For instance, 
in ^ the collective label j* stands for ji, j2, js and j4, where each strand gets its own 
spin label. The blank and primed circles convey whether it is the group element g or 
g' that is taken in the given representation. The horizontal dotted line represents the 
triangulation tetrahedron dual to the dual edge to which the given group element is 
associated. 

The first equality in Q follows directly from the normalization of the Haar group 
measure, its invariance under translations and the multiplicative property of represen- 
tation matrix elements. In this context, group integration is also known to produce 
a projection operator onto the space of intertwiners among the representations given 
by the four strands. The last equality in ([6| illustrates this identity by expanding this 
projector in a basis of normalized intertwiners; the bracketed spin network provides 
correct normalization in the denominator of the expression. The summation over the 
new intertwiner basis labels i and i' make up the sums over dual edge labels, part of 
the sum over spin foams in Q . Performing the same group integration and intertwiner 
expansion over all dual edges of the triangulation, we can read off the edge and vertex 
ampHtudes of equation Q . 

Thus, for discrete BF theory, writing all tensor contractions in terms of spin 
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networks we find these amplitudes to be 




The topology of the contraction graph corresponding to above follows directly 
from the adjacency structure of the dual 2-complex. We shall refer to this graph as 
the pent graph; it will appear in the vertex amplitude definition of each spin foam 
model discussed later in this section. Both the edge and vertex amplitudes, and 
Ay, appear with full spin labelling. For conciseness, most of the spin labels will be 
suppressed or represented schematically, as in equation ([6]), in the rest of the paper. 

Starting from this basic setup, new models may be obtained by modifying the 
partition function directly, by changing amplitudes at the level of equation Q, or at 
an intermediate level, by modifying the integrands in equation ^ . We present the new 
models following the last approach. Note that the subsequent discussion emphasizes 
not the derivation of these models from first principles, but their formulation in a 
unified framework suitable for computation. 



2.1 Gravity, Barrett-Crane and new models 

The Barrett-Crane (BC) model starts with the quantized BF theory path integral Q 
and imposes restrictions on the spin foam summation in equation Q. These re- 
strictions heuristically correspond to imposing the simplicity constraints on the B 



field 10 . The restriction is twofold. First, the Spin(4) representations are restricted 
to balanced ones j = where = j+ = j. Second, the intertwiner summation 

and edge weights of equation (|6| are modified such that the z-sums contain only a 
single term corresponding to the so-called BC 4-valent intertwiner. 



The BC model amplitudes are given in section 2.2.1 The evaluation of this vertex 



amplitude is discussed in several papers [8|[9 16 ,25 , where variations on the face and 
edge amplitudes have also been considered. 

Recently, shortcomings of the BC model have been identified by several authors. 
Specifically, while this vertex amplitude correctly reproduced the asymptotic behavior 
of some graviton propagator matrix elements, it does not do so for all of them [3j|4j[34]. 
Modified spin foam models, referred to here as new models, have been subsequently 
proposed with the hope of overcoming these difficulties. The model proposed by Engle, 
Pereira, and Rovelli (EPR) [20 21 and by Livine and Speziale [27][28| had the common 
motivation of identifying its boundary state space with the space of spin network states 
of loop quantum gravity. The model proposed by Freidel and Krasnov (FK) [22] was 
derived in a similar fashion, but made different choices while imposing the simplicity 
constraints. As a result, the FK model's boundary state space is different from that 
of the EPR one. More recently, Conrady and Freidel have discussed in more detail 
the boundary state space of the FK model [Tsl. 
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2.2 Model framework 



The BC, EPR, and FK models may be presented within the same framework, fol- 
lowing [22) . We briefly present this framework and how each model is realized in 
it. 

The first step, compared to BF theory, as for the BC model, is to restrict the 
Spin(4) representations to balanced ones, j = 

Consider a single strand from the double integral in equation It depicts the 
product of two linear operators, corresponding to group elements g and g' , in the 
Spin(4) irrep j. One could always insert the identity operator between g and g' without 
changing anything. On the other hand, inserting a different linear operator in the 
same place will produce different results. Keeping with the goal of identifying the i- 
intertwiners in ^ with the SU (2) intertwiners in loop quantum gravity spin networks, 
instead of inserting an arbitrary linear operator, we insert only an arbitrary SU{2) 
intertwiner. The difference between the spin foam models described in the following 
sections comes down to the choice of this insertion. An SU (2) intertwiner is inserted 
in the following sense. 

Because of the decomposition Spin(4) = SU (2) x SU (2) , a balanced irrep of Spin(4) 
can be written as the tensor product of two copies of an SU (2) irrep, j = (j, j) = 
Seen as a tensor product of two SU{2) reps, using Clebsch-Gordan rules, it can be 
decomposed into a direct sum of SU{2) representations 0, 2, . . . , 2j. This decomposi- 
tion is not unique, since one is possible for each diagonal S'?7(2)-subgroup of Spin(4). 
However, each such choice is equivalent because of the group integrals surrounding 
the operator insertion [cf. ([6])]; ultimately, spin foam amplitudes are independent of 
the choice. By inserting an SU(2) intertwiner, j ® j ^ j ® J, between g and g' , 
we essentially insert a linear combination of projections onto the SU(2) irreducible 
invariant subspaces: 



^ O 2fc 



2k 



o 



2k 



, (8) 



where the Cjk are weight factors parametrizing the choice of intertwiner and the 
trivalent vertex corresponds to the Clebsch-Gordan tensor. The extra coefficients are 
necessary, in our choice of normalization, for the Clebsch-Gordan decomposition of 
the j ® j representation. 

Each strand in the above diagram carries a spin label. For compactness of pre- 
sentation, some labels are listed separated by commas, next to a network (instead of 
being directly attached to the corresponding strand) or even omitted. Some groups of 
spins are also represented by a collective label like j* . No ambiguity should arise, as 
the spin labelling is essentially unique, given the equality relations. This convention 
is used throughout the rest of the paper. 

The overall spin network conventions and normalizations used in this paper are 



those of 13 24 . Their precise relation to S'[/(2)-tcnsor contractions is presented in 
detail in the Appendix of 14 . Note that the intermediate spiij^ is always even and 



^The notation of this paper differs from reference |22|, as Freidel and Krasnov use half-integral 
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so can be written as 2k for integer k. 

Consider just a single Clebsch-Gordan projector insertion for each of the strands 
in (|6]) , as show in ([s]) , and concentrate only on the part of that equation below the 
dotted line. The same calculation will have to be done above the dotted line, that is 
on, both sides of each tetrahedron in the triangulation. The Spin(4) group integral 
can be written as an integral over two copies of SU(2), where g — {g^,g'^)- The 
matrix elements of g in a balanced representation j will be an outer product of matrix 
elements of .g~ and g+ in representation j (circles with — and +, respectively): 

2k* 



2k 



I 



dg dg^ 





(9) 



Summations over the intertwiners again follow directly from the property that 
group integration is equivalent to projection onto the space of intertwiners between 
the four j representations (j -spins). The extra summation over the i intertwiners can 
be inserted because the Clebsch-Gordan projectors map each pair of intertwiners 
into the subspace of intertwiners between the four 2k representations (k-spins). These 
intertwiners can be conveniently parametrized, as depicted, by an even integer 2i (i- 
spins^ 

The open spin networks at the bottom of the right hand side of ^ join with 
other similar spin networks and form left (— ) and right (+) pent networks, which 
will contribute to the corresponding vertex amplitude. These are spin networks with 
topology shown in ( |33[ ), obtained by substitution of z~ and i"*" intertwiners into the 
pent graph of equation The open spin network at the top of the same expression 
joins with its mirror image above the dotted line and contributes to the corresponding 
edge amplitude. With the exception of the tripetal spin networl^locaied at the center 
of the diagram on the right of all spin networks appearing so far have known 
evaluations. They have come up in the evaluation of the BC vertex amplitude and 
have been explicitly computed using recoupling techniques from 13 24 . The tripetal 
spin network will be evaluated in section |4j 



spins, while we use integral twice-spins to label SU{2) irreps. In the present notation, their j becomes 
j/2, but their k does coincide numerically with our k. Also, their / corresponds to i introduced in 
equation |9]l. 

*It should be noted that reference [20| uses half-integral spins, while we use integral twice-spins 
to label SU (2) irreps. However, Engle, Pereira and Rovelli's definitions for i- and j'-spins numerically 
coincide with ours. 

^This spin network was first introduced in equation (5) of 20 and is also referred to as a fusion 
coefficient. The name tripetal is suggestive of the topology of the network, in which the spins i"*", i~ 
and 2i correspond to the edges of the three petals. 
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To completely define each of the three models under consideration, it remains only 
to specify the dual 2-skeleton spin foam labelling and the weight factors Cjk in 
The face, edge and vertex amplitudes are then specified by the preceding construction 
(see section |4] for an important caveat). Most generally, in this framework, the spin 
foams summed over in the partition function Q assign a j-spin to each dual face, an 
«-spin to each dual edge, and a fc-spin to each dual edge-dual face pair. However, the 
number of labels may be reduced in particular models. At present, only the FK model 
uses all label types. 

As mentioned earlier, while this presentation is convenient for computational pur- 
poses, it hides some of the motivation from the derivation of these models. More 
physical insight for each model can be found in the original references. 

2.2.1 BC model 

In the Barrett-Crane model, the faces of the dual 2-complex are labelled by j-spins. 
The choice of intertwiner insertion weights are Cjk — Q{k^Q), and Cjo = (— )"' (j + l); 
which is the value of the j-loop spin network. Each dual edge is shared by 4 dual faces, 
while each dual vertex is shared by 10 dual faces. The preceding construction specifies 
the following dual face, edge and vertex amplitudes: 



^/(^) = (i/ + l) , MF) = l. and A,{F)= W \ X / ' (1^^) 




where 





(10b) 



Here the spin arguments are determined by the dual faces sharing the given cell of the 
2-skeleton. Specifically, the vertex amplitude depends on 10 spins, hence its name, 
the BC lOj-symbol. It is important to note that different edge and face amplitudes 
have been proposed for the BC model as well [9 19 32 . 



2.2.2 EPR model 

In the Engle-Pereira-RovcUi model, the dual faces are labelled by j'-spins, and dual 
edges are labelled by i-spins [cf. ([9])]. The weights are Cjj — 1 and Cjk — for k ^ j. 
Each dual edge is shared by 4 dual faces, while each dual vertex is shared by 10 dual 
faces and 5 dual edges. The preceding construction specifies the following dual face. 
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edge and vertex amplitudes: 

Af{F) = {jf + ir, A,{F) = 



O 2i, 



^3 ^^2j*,2le 



(11a) 



and 




(lib) 



where 



O.- Oi+ 




(11c) 



Here the spin arguments are determined by the dual faces and edges sharing the given 
cell of the 2-skeleton. Specifically, the vertex amplitude depends on the 10 j-spins 
from the dual faces sharing it, as well as the 5 i-spins from the incident dual edges, 
hence it may be called the EPR 15j-symbol. The same vertex amplitude was derived 
and 



20 



28 , although the former reference was not specific about face and edge 



amplitudes. 



2.2.3 FK model 



In the Freidel-Krasnov model, the dual faces are again labelled by spins, denoted j, 
and dual edges are also labelled by intertwiners, denoted i, and finally each dual edge- 
face pair contributes an independent spin, denoted k [cf. ([9|]. The weight factor is 
more complex than for the other two models and is given by 



C. 



jk 



(j-fc)!(j + fc + l)!- 



(12) 



Each dual edge is shared by 4 dual faces, while each dual vertex is shared by 10 dual 
faces and 5 dual edges, as well as 20 individual dual edge-face pairs. The preceding 
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construction specifies the following dual face, edge and vertex amplitudes 

O 2Ze 



AfiF) = Uf + ir, A,{F) 



and 




(13b) 



where 




(13c) 

Here the spin arguments are determined by the dual faces and edges sharing the given 
cell of the 2-skeleton. Specifically, the vertex amplitude depends on the 10 j-spins 
from the dual faces sharing it, as well as the 5 i-spins from the incident dual edges, 
and on the 20 fc-spins from the dual edge- face pairs sharing it. Thus it may be called 
the FK 35j-symbol. Setting all fc-spins, and necessarily all i-spins, to 0, this model 
exactly reproduces the BC spin foam amplitudes. Also, setting all fc-spins equal to the 
corresponding j-spins exactly reproduces the EPR vertex amplitude Ay{F). However, 
in that case, the EPR edge amplitude Ae{F) is reproduced with the extra factor 

n/3e[(j/ + l)!]V(2j/ + l)!. 

Here is a brief summary of the spin labelling conventions: The BC model assigns 
integer labels only to dual faces {j -spins). The EPR model also assigns integer labels 
to dual edges {i-spins). The FK model additionally assigns integers to each dual edge- 
dual face pair (k-spins). The full labelling is illustrated explicitly for a single pent 
graph in figure [T| A detailed explanation of label notation follows at the end of the 
next section. 



2.3 Observables and boundary states 

If the underlying triangulated manifold is closed, then corresponding spin foams are 
also said to be closed. Similarly, if the underlying manifold has a boundary, the spin 
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foams are said to be open and also have a boundary. Any open spin foam Fq can 
be decomposed into Fq = F Li Fb, where Fb labels only cells dual to the boundary, 
while F labels only cells dual to triangulation simplices in the interior. For an open 
spin foam Fq, its amplitude may be naturally generalized to 

AiFo)^ A{F,FB)^iFB), (14) 

where the bulk amplitude A{F,Fb) is the usual amplitude defined according to 
and ^! is referred to as the boundary state, which may be fixed separately from the 
bulk amplitude. The partition function in the presence of a boundary is then written 
as 

Z^= MF,FB)^iFB). (15) 

F,Fb 

Observables are functions 0(F) defined on spin foams. Their expectation values, in 
the absence and in the presence of a boundary, are respectively 

(O) = ^Y.OiF)AiF) and (O)* = ^ ^ 0(F, Fb)A(F, Fb)vI/(Fb). (16) 

F * F,Fb 

As an illustration, an open spin foam model with a boundary state may arise 
if we split a closed spin foam model in two parts and average over one of them. 
Suppose a close triangulated manifold can be decomposed into two bulk pieces and 
the codimension-1 boundary between them. Any closed spin foam Fc can then be 
decomposed as Fc = F U Fb U F', where Fb corresponds to the boundary, F to the 
interior of the piece we are interested in and F' to the interior of the other piece. The 
partition function may be rewritten as follows: 

Z = J2 MF, Fb)A{Fb)A{F', Fb) 

Fc 

= J2 MF, Fb)A{Fb) J2 MF', Fb) = ^ A(F, Fb)-^{Fb). (17) 

F,Fb F' F,Fb 

where "^(Fb) has been defined by averaging over all spin foams F'. This example is 
very similar to the separation of a large system into a subsystem and the environment 
in quantum statistical mechanics. 

The simplest example of a triangulation with boundary is a single 4-simplex, with 
the five tetrahedra forming its boundary. The 2-complex dual to the interior consists 
of a single dual vertex, corresponding to the 4-simplex itself. The dual 2-complex of 
the boundary consists of five dual edges, dual to the tetrahedra, and of ten dual faces, 
dual to the triangular faces of the tetrahedra. The problems described in sections |3.2| 



and 3.3 have previously only been considered for a single 4-simplex. This paper 
restricts attention to the same case. 

The algorithms that will be described in section |5] are apphcable only to a restricted 
class of states, factored states. Such a state must factor in a specific way with respect 
to the spins it depends on. The various spin labels of the dual complex of the 4-simplex 
and the corresponding notation are summarized on the pent graph of figure [l] The 
vertices of the pent graph correspond to the five boundary tetrahedra of the 4-simplex, 
while the ten edges connecting them correspond to its triangles. This graph is labelled 
by 35 spins, i^, jj,,ei ^^nd k^^. The e subscript numbers the vertices of the pent graph; 
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Figure 1: The pent graph, summarizing the indexing scheme for i-, j-, and A:-spins. 



it is always taken mod 5. The spin jc,e labels the graph edge joining vertices e and 
e + c. The superscript x stands for either p or q; ^ labels the vertex-edge pair e and 
(c, e), while k^^ labels the pair e + c and (c, e). Again, all vertex indices are taken 
mod 5. 

The class of factored states is somewhat different for each model. However, it 
contains at least all of the following: 

= n ^cAjc.e) n n ^-j-Akij, m 

c,e e x,f,e 

where the products range over all i-, j-, and fc-spins. Spins not part of a particular 
model may be dropped from the product. The ips are arbitrary functions with finite 
support (more on that below). For each model, the class of factored states is enlarged, 
as factors of \1/(Fb) may be allowed to depend on specific clusters of spins, instead of 
only individual ones. The details will be elaborated in section |5] 

Nearly all previous work on the problems described in sections |3.2| and |3.3| has 
considered only factored boundary states. While this class of states is restrictive, its 



limitations may be overcome. Note that the expectation value (0)>p in equation (16 1 
is equal to the ratio of two quantities that are both linear in the boundary state 
The numerical algorithm computes this numerator and denominator separately. Since 
any boundary state \[' can be approximated by linear combinations of factored states 
with finite support, so can (O)* be approximated for any boundary state However, 
the efficiency of such an approach is yet to be analized. The condition of finite support 
for the factors i^U) is crucial for the sums defining (O)* to be finite. 



3 Observables and asymptotics 

One of the motivations for constructing new vertex amplitudes is the recently discov- 
ered inadequacy of the BC model in reproducing semiclassical graviton propagator 
behavior in the large spin limit. Some of the propagator matrix elements show the 
expected behavior, while others do not [3]|^[34]. Thus, it is important to identify 
where the new vertex models differ from the BC one and whether they have better 
semiclassical limits. 

The comparison should ultimately be done at the level of physical observables 
computed within each model. An important class of observables, already mentioned 
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above, are matrix elements of the graviton propagator. They will be discussed below in 



section 3.3 Another class of observables associated with propagation of semiclassical 
wave packets is addressed in section |3.2[ Finally, a possibly less physically meaningful 
but technically simpler comparison can be made at the level of vertex amplitudes. 
It too can reveal important information about the behavior of the models and is 
described first below. 

3.1 Comparison of asymptotics 

One complication is the difference in the spin argument structures: the BC vertex 
has 10 spin arguments, the EPR vertex has 15 spin arguments, while the FK vertex 
has a total of 35 independent spin arguments. This complication may be overcome by 
fixing the 10 common j-spin arguments and maximizing the vertex amplitude over the 
remaining spins. This effective vertex amplitude can be substituted into the partition 
function Q where the summation is then performed over spin foams which only 
assign j-spin labels to the dual 2-complex. This simplification allows the comparison 
of amplitudes for individual spin foams. 

It is important to note that the amplitudes in Q have contributions from faces 
and edges as well as vertices. The face amplitudes are the same for all models and 
are easily factored out. The edge amplitudes, on the other hand, also differ from 
model to model and thus must be included in the amplitude comparison. To make 
the comparison on a vertex by vertex basis, the edge amplitudes are split between the 
vertices they connect as follows: 



(j) ~maxA,(j,fc,j) m \A,{j,k,i)\. (19) 

' y eD-u 

For simplicity, we consider only spin inputs where each of the jc,e, ie and fc^g sets of 
spins have equal values, respectively denoted by j, i, and k. Our assumption is that 
vertex amplitudes for these spin inputs behave generically. Small scale numerical tests 



support this assumption. Otherwise, maximizing the expression in ( 19 1 over a larger 
i, fc-parameter space quickly becomes impractical. 

3.2 Semiclassical wave packets 



The problem presented in this section was introduced in 29 . Consider a single 4- 
simplex. As shown in the preceding section, it is described by a spin foam with a single 
dual vertex and i-, j-, and fc-spins labelling cells dual to its boundary. An arbitrary 
functional ^!{Fb) depending on these boundary spins, in general, corresponds to a 
statistical quantum state, that is, a density matrix. 

This is analogous to the single point particle, where an arbitrary density ma- 
trix p{xf,Xi) — {xf,tf\p\xi,ti) can be described in terms of its matrix elements be- 
tween eigenstates of the Heisenberg position operator at different timesj^ x{ti)\xi, ti) — 
Xi\xi,ti) and x{tf)\xf,tf) — Xf\xf,tf). The density matrix is pure only if it can be 
factored, p{xf,Xi) = ^{xf^tf)*^j{xi,ti), where ^{x,t) denotes the time evolution of a 
given wave function. 



®In this representation, the functional p{xf,Xi) is not necessarily symmetric, p{xi,Xf) 7^ 

p{Xf,Xi)*. 
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Similarly, we can split the boundary of the 4-simplex into two piece^ the initial 
(— ) and the final (+). Then, for a pure boundary state, we should be able to write 

nFB) = ^+iF+r^-iF^), (20) 

where respectively depend only on spins labelling the dual complex of the corre- 
sponding piece of the boundary. 

The relationship of the two boundary state factors '^±{Fg) is constrained in two 
ways. On the one hand (in the limit of h ^ 0), the amplitude should be peaked on 
those geometries that correspond to the boundary of a classical 4-geometry satisfying 
Einstein's equations. On the other hand, should be a time-evolved, "future" 
version of the "past" which can be expressed as 

^+{F+) = J2 A{F,Fb)^-{F^), (21) 

F 

where the summation over the boundary spin foams keeps Ft fixed and varies Fj7 . 



Reference 29 has proposed an expression for ^'(Fb), in the context of the EPR 
model, which should reproduce a fiat regular 4-simplex. This state has gaussian 
dependence on individual spins and hence is factorable in a convenient way. The 



problem is then to compute ^'_|_(F^) both from (20 1 and from (21), and to compare 
the two. Agreement is interpreted as evidence of a correct semiclassical limit for the 
EPR model. 

A concrete expression for the proposed "^{Fb) is 

■^{FB)=N\{^,^,{j,^,)\[i^SeAjcA). with (22) 

c.e e 

V'c,e(jc,e)=e-*(>--^«)'+^^^-, (23) 



i^SeAjcA) = \ L,,- wr^. TT. .e-w(--)^+»t-, (24) 

y 9{2le,2ji^e, 2J2,e)t'(2^e, 2ji^e-l, 2:;2,e-2) 

where TV is a normalization factor, jo determines the size of the regular 4-simplex and 
cosO = —1/4. The parameter r controls the size of quantum fiuctuations about the 
classical values of j. 

The wave packet propagation geometry given in f29l fixes r = in the state ( 23 1 , 



freezing all j-spins to the background value jq. Effectively, only the dependence of 
'^{Fb) on the i-spins was considered. A single vertex of the 4-simplex is labelled 
as "past" , while the remaining four as "future" . The four "future" vertices form a 
tetrahedron, whose dual is labelled by an i-spin. This labelled dual edge constitutes 
Fg, while the remaining four dual edges labelled by i-spins constitute Fg. This 
propagation geometry will be referred to as EPR 4-1 propagation. 

An immediate generalization, feasible with the algorithm described in section [5] is 
to relax the r = limitation. The choice of t should be consistent with the parameters 
used in the graviton propagator calculations. Thus, following ^l7j, we let the wave 

'^Technically speaking, this decomposition is unique only in Lorentzian signature. In Rieman- 
nian signature, different choices of the decomposition should correspond to different possible Wicls 
rotations. 

*It should be noted that reference [l7| uses half-integral spins, while we use integral twice-spins 
to label SU{2) irreps. A j label from Christensen, Livine and Speziale corresponds numerically to 
j/2 in current notation. 
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packet width depend on the background spin, 

T = 4jo/a, (25) 

with a is a positive parameter. 



It is important to note that the boundary state proposed above, in equation (22 1, 
is not the best possible candidate to generahze the calculations of [29]. Its major 
advantage is that it belongs to the class of factored states. The major disadvantage 
of factored states is that they neglect possible spin-spin correlations. In fact, a more 
realistic proposal for a boundary state where both i- and j'-spins are allowed to vary 



and which includes such correlations was given by Rovelli and Speziale 35 . Unfortu- 
nately, Rovelli and Speziale 's boundary state does not factor nicely and would require 
more sophisticated techniques to be used efficiently. As such, the boundary state pro- 
posed above should be seen more as a test of the algorithms presented in section [5] and 
an attempt to explore the qualitative effects of introducing j-dependent (and below 
fc-dependent) boundary states. 

Now, a uniform methodology should be constructed for each of the three models. 
As only j-spins are common among the models, we propose the following wave packet 
propagation geometry. One possibility is to propagate wave packets from nine of the 
j-spins to the remaining one. This configuration corresponds to fixing a single triangle 
(defined by three vertices of a 4-simplex) in the "future" , while relegating the other 
nine triangles (containing at least one of the two remaining 4-simplex vertices) to the 
"past" . Thus, the single j'-labelled face dual to the "future" triangle will constitute 
Fg , while the rest of the boundary spin foam will constitute , including all i- and 
fc-spins, if any. This propagation geometry will be referred to as 9-1 propagation. 

Another alternative is to assign a vertex of the 4-simplex to the "future" , together 
with the six triangles sharing sharing it. The remaining four triangles are relegated to 
the "past" . Thus, F^ consists of the six j-labelled faces dual to the "future" triangles, 
with the rest of the boundary spin foam constituting Fg . This propagation geometry 
will be referred to as 4-6 propagation. There are numerous other possibilities. However, 
the two described above are sufficient to illustrate an application of the numerical 
algorithms and to show the qualitative behavior to be expected from propagated wave 
packets. 



The boundary state (|22| is valid only for the EPR model. For the BC model, we 



simply drop the tpe factors: 



And for the FK model we must add extra V'c e factors for each /c-spin: 

^(Fb) - Nlli,,,,{j,,,)l[M^e, {Jce}) n ^cAK,e,Jc,e)- (27) 
c,e e x,c.e 

Because the fc-spins are closely geometrically associated with j-spins, we use the same 
gaussian state parameters: 



U-ky.u + k + iy.- 
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The square root factor includes the FK model edge normalization, as does (24) for 
the EPR model. 



3.3 Graviton propagator 

The graviton propagator is well defined in the perturbative quantization of gravity. 
It is computed as the 2-point function G^vpa{x,y) — {Q\h^j_y{x)hpa{y)\Q) , where |0) 
is the Minkowski vacuum, and h^u^x) is the metric perturbation. General relativity 
requires that, in harmonic gauge |37], the decay rate of the 2-point function, for 
large separation between points x and y, is the same as for the Newtonian force 
of gravitational attraction: inverse distance squared. The framework for computing 
the equivalent of the graviton propagator in the spin foam formalism was elaborated 
in [T2[[26}|30}[34| . The quantum area spectrum is A = ^p(j + 1), with j a dual face spin 
foam label and ip the Plank length. Dimensional arguments then give the expected 
decay of the propagator as 0(l/j), with j being the typical size for the chosen spin 
foam boundary state. 

The expected asymptotic behavior of the graviton propagator has been checked 



for the BC model both analytically and numerically 17 26 34 . Unfortunately, the 
expected behavior was only reproduced for certain tensor components of G'^^po-(x, y), 
but not for others [3j[4] . This negative result has prompted the introduction of EPR 
and FK spin foam models as alternatives to the BC model. The challenge is to 
compute the graviton propagator for the new models and check that it has the expected 
asymptotic behavior. 



Following 17 , we show the computational set up for the BC model and then 
extend it to other models. Consider again a single 4-simplex with boundary and the 
corresponding spin foam. We associate the area A — ^p(j + 1) to each triangle, 
depending on the j-spin labelling its dual. The goal is to compute the correlation 



between observables depending on the triangle areas [cf. ( 16 1] 



Wce,c'e' ^^y\ A{F,FB)h,,{FB)h,,,,{FB)^{,FB), (30) 
^* F.F, 

where ce and c'e' index the specific jc,e and jc' ,e' spins taking part in the correlation. 
Again following 17 , the boundary stat^ is a semiclassical gaussian state peaked 



around a flat 4-simplex, whose scale is set by jq: 

^{Fb) = []e-3^(^-="^°^'+*®^-% (31) 

c.e 

where cos0 = —1/4, and jo sets the scale for the background geometry, as in the 
previous section. The observables measure the fluctuation of areas squared: 

h (F \- O'c.e + 1)' - (jO + 1)' 

hUFB) - ^r—^, . (32) 

Note that the product "^'{Fb) — hce{FB)'i^ {Fb) has exactly the same factorizability 
properties as ^>(Fb). This property allows both the numerator and denominator 



in (30 1 to be computed on the same footing. 



®We incorporate the "measure" discussed in [17| into the boundary state and pick the trivial case 
k = 0. 
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Again, an important task here is the generaHzation of this calculation to the EPR 
and FK models. This generalization essentially requires the specification of a boundary 
state that describes a semiclassical state peaked around the flat regular 4-simplex. 
Since this is the same requirement used in picking out the boundary states in the 
section on wave packet propagation, simply choose the same ones. That is, the BC, 



EPR and FK boundary states are specified, respectively, by equations ( 26 1 , ( 22 ) , and 



(271 



4 Spin network evaluations 

The second group integration identity in (|6| requires a choice of basis in the space of 
intertwiners between four Spin(4) representations. This choice is arbitrary; however, 
some choices are more convenient than others. For example, the normalization factor 
in ^ is simplest when the i'-basis is the same as the i-basis. On the other hand, 



the choice of intertwiner basis in the EPR model, equation (11c I, is made such that 
when the intertwiner networks are substituted into the vertex amplitude pent graph, 
equation (|7]), the amplitude is resolved as a sum over 15j-symbols with the following 
topology: 




(33) 



A similar choice is made for the BC and FK models, equations ( 10b[ ) and (13c). This 



topology is required for the numerical algorithm described in section [sj 

Unfortunately, the requirements of simple edge normalization factors and the above 
topology requirement for each vertex amplitude are not always compatible. For ex- 
ample, they are not compatible for the minimal triangulation of the 4-sphere. With 
the current formulation of the vertex amplitude evaluation algorithm, preference must 
be given to the topology requirement. A similar issue, referred to as "edge splitting," 
was encountered for spin foams on a cubic lattice in [14| . 

Throughout this paper, we have assumed that the topology and simpHcity of dual 
edge normalization requirements can be simultaneously satisfied. This assumption 
is justified in the case of a single 4-simplex, and other simple arrangements of a 
small number of 4-simplices. If this assumption is not justified, then the dual edge 
amplitudes given in the previous section will have to be modified, with the important 
exception of the BC model. The edge normalization requirement for the BC model is 
trivial. 

4.1 Tripetal network evaluation 

The tripetal spin network, defined in equation ([9]), is evaluated as follows. It is first 
written as a contraction of two trivalent networks, along the strands labelled , i+ 
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E 





j j f 




■/ / j' 




'2k' 


2k n 


i~ n 2k 




i+ n 2k' 






i+ 2i 




6{i-,n,2k) 


e{i+,n,2k') 


9{i- 


,i+,2i) 




(37) 



The first equality recouples the crossing strands through the auxiliary spin n. Other 
steps correspond to collapsing triangles to 3-valent vertices. The square brackets 
denote the evaluation of a tetrahedral (tet) spin network, while 6 and A denote the 
evaluations of the theta-like and loop spin networks seen in (|8| and elsewhere. For 



full details, see 24 and the Appendix of 25 



After this simplification, the tripetal network is proportional to the theta network, 
where we write the left copy of the above coefficient as Pi~i+ and the right copy as 
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Q^ 




2i 




(38) 



The network on the right hand side of the above equation is equal to the theta network 
up to sign, which is specified in (39a I. Both Pi-i+ and Qi-i+ depend on many spins. 



The displayed indices are those that will be important in section |5] 

To obtain the final formulas for each vertex e of the pent graph, we make appro- 
priate substitutions into the above expression, from each half of the tripetal network. 
We replace i and by ie and if respectively. In the P coefficient, the spin k becomes 
fc^ e; while k' becomes fcj^. At the same time, the spin j becomes ji^e, while j' be- 
comes j2,e- Similarly, in Q, the spin k becomes ^2 6-2' ^' becomes k^ j becomes 
j2,e-2 and j' becomes ji.e-i- The final formula for the tripetal network is 



R': 



(39a) 



2L 



Jl,e 



J2,e 



J2,e 
i + 



J2,e 



Jl,e 
2^2.6 



&{np,ji^e,j2,e) e{ie,np,2kU e{it,np,2k^J 



(39b) 



hU2 



-2+Jl,e-l-" 



■^'^2,e-2 



Oinqj2..e-2jl,e-l) 



j2,e-2 j2,e-2 Jl,e-1 

e{ie ,nq,2kl^_^) 



Jl,e-1 Jl,e-1 j2,e-2 



n„ 



2h'' 



{it ,nq, 2k1 



(39c) 



Notation for the j- and fc-spins is explained in the next section. 

The same spin network has also been evaluated independently by Alesci et al [2], 
however their notation differs from ours in normalization. 



5 Numerical algorithms 

The dual face and edge amplitudes, Af and A^, are trivial to compute. The difficulty 
lies in evaluating the dual vertex amplitude Ay, which is where we will concentrate. 
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All algorithms described in this section are extensions of the original CE algorithm 
for the BC model [l6] and all fall into the same product-trace pattern: 

A.({je.e,ie,fc,%}) = ^ tr[M4M3M2Mi Afo] , (40) 

,m+ 

where {—)^ is an overall sign factor, (f) depends only on to* and are matrices of 
compatible dimensions, collectively depending on all the spins. Each of these elements 
may be specified separately in any incarnation of this algorithm. In all cases presented 
below, we have 

(j) ^ (-)5(™"-™+)(to- + l)(m+ + 1). (41) 

However, the Mg matrices will be redefined for each variation of the algorithm. The 
notation for various boundary spins is summarized with the pent graph in figure [T| 

The reason for the similarity is that all the vertex amplitudes discussed in section [2] 
are structurally similar to the vertex amplitude ([t]) for discrete BF theory, in which the 
intertwiners can be resolved to produce SU{2) spin networks of the topology shown at 
the beginning of scction|4] The topology of the pent graph, as described in flBl, is what 
gives rise to the product trace structure pO| ) . The differences between the Mg matrices 
for each model capture the differences between the Spin(4) intertwiners used to define 



the corresponding vertex amplitudes in equations (10b I, (11c I, and (13c I. When 
contracting with a boundary state, the structure of the matrices also captures the 
extra summations with respect to the boundary spins. 

5.1 Run time complexity 

The discussion below assumes familiarity with run time complexity estimates of stan- 
dard numerical linear algebra operations. 

The run time complexity of a generalized CE algorithm may be estimated as fol- 



lows. Suppose that the spin arguments to Ay in (40 1 are of average magnitude j. 
Then, generally, the dimensions of the matrices scale as a power of j; say, each 
matrix is 0{j'^) x 0{f^), for some integer d. The run time will be dominated by fill- 
ing the Mf, matrices and by the product-trace operation. This exponent d must be 
determined separately, depending on the model and on the inclusion of a boundary 
state. 

The product-trace may be implemented as follows: each of the 0{j'^) standard basis 
vectors is subjected to matrix-vector multiplies by the and appropriate elements 
of the result vectors are accumulated into the trace. If the Me are dense, then the 
cost of a matrix- vector multiply is 0{j^'^). However, we shall see below that this 
complexity may be reduced by decomposing each into sparse factors. Hence, we 
will parametrize the matrix-vector multiply complexity as 0{j^), with D no greater 
than 2d, and the product-trace complexity as 0{j'''~^^). 

If each matrix element is computed in at most 0{i^) time, the upper bound on 
the time needed to fill an 0{j'^) x 0{j'^) matrix Me is 0{j^'^^^). Sparse factorization 
improves this estimate as well, which we will parametrize as 0{j^^f), where F also 
does not exceed 2d. In all cases we have examined, d + D > F + f, which implies that 
the product-trace operation dominates matrix filling in run time for large spins. More 
detailed discussions of possible optimizations for matrix filling can be found in ^25]. 
Below, we will give the best known value of / for each algorithm. In particular, if a 
matrix element depends on a tet, the standard Wigner-Racah formula requires 0(j) 
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operations to compute it (/ = 1). On the other hand, using recurrence relations [36] 
or hashing techniques the number of operations can be reduced to 0(1) (/ = 0). 

Finally, the outer sums in (40) also span ranges of size 0{j). Therefore, the 
run time complexity of a generalized CE algorithm may be expressed as 0(j^~'"'*+^). 

The matrix elements of the (computed for each case in the following sections), 
will contain spin network evaluations that require certain inequality and parity con- 
straints on their arguments. Solving these constraints yields precise matrix dimensions 
and bounds for any intermediate summations. The details are described at the end of 
section I5T2I It can be shown that all matrix dimensions as well as intermediate sum- 
mation bounds are finite. However, for that to be true in the presence of a boundary 
state, it is crucial that each factor of the boundary state (18 1 has finite support. 



5.2 Algorithms for fixed boundary spins 
5.2.1 BC vertex 

For the BC model, the vertex amplitude is only a function of the j-spins. As a slight 
abuse of notation, we will use the symbols ig as indices (also referred to as spins, 
and directly analogous to the if indices introduced for the other models) of the Mg 
matrices: 



(Me) 



{ie + 1) 


ie j2.e m 




ie J2,e 111+ 




ie+1 32,e-l jl,e 




ie+1 j2,e-l jl,e 





Q{j2,e-l,ie+l,jl,e) ^'(j2,e , «e , jl,e) 6'( j2,e , ^e+l , ™ ) 0{j2,e, ie+1 , m+) ' 



(42) 



The sign factor from (40 1 is given by S" = J2cejc,e- The ranges of the ig and m spins 
are specified by triangle inequalities and parity constraints satisfied by various spins. 
For a detailed derivation and for notation, see the original reference 16 , and also |25] 
and its Appendix]^ 

The structure of the matrices will become increasingly important and will 
grow in sophistication in the algorithms presented below. Hence, it is convenient to 
introduce a graphical notation to represent this structure. In this simplest case we 
have: 



(43) 



Each strand represents an index. The incoming and outgoing strands correspond to 
the ie and ie+i indices of Mg and are labelled as such. The product-trace operation 
in (40) is effected by concatenating appropriately labelled strands. Further features 



of this graphical notation will be elaborated as they are introduced. 

Each matrix Ale is dense and of size 0{j) x 0{j). According to the discussion at 
the beginning of this section, we have d = 1, _D = 2, / = 0, F = 2, and 2 + d + D = 
5. Hence, we recover the well known O(j^) run time complexity of the original CE 
algorithm. 



5.2.2 New vertices 

For the EPR model, the amplitude is a function of both i- and j-spins, while the FK 
model is also a function of fc-spins. Here, we give the explicit FK formula, with the 

'^"Note however, that these references use half-integral spins, while the present paper uses integer 
twice-spins. 
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EPR version obtained by setting fc^^ 



The matrix elements of are 



+ l*e+l *e i 



P'L .+ Ar^_ ivv ^ 

^^^^ »e ^. d{le,lt,2l,) 



(44) 



where 



+ 1' 



J2,e-1 
^=^,i2,e) 



(«?,Ji, 



,j2,e)6'(«?,j2,e-2, jl,e-l) ' 



(45) 
(46) 



while P and Q are given by equations ( 39b ) and ( 39c ) . A detailed discussion of the 



bounds on the various spins follows below. The T± and N terms are inherited from 



the CE algorithm; NT^T+ is a particular factorization of the right hand side of (42 1 



The new factors of P and Q come from the evaluation of the tripetal networks in 

,«Jj_i). In 



equations (11c I and (13c 



The matrix elements of Mg are indexed by the pairs {i^ and {ig_^_i, tg_^_i 
graphical notation, Mg has the following structured factorization. 



■■e+l 



-e+1 




(47) 



where P*^ stands for the entire bracketed term in (44 1. Unmarked vertices in the above 



diagram essentially correspond to Kronecker deltas. The notation is saying that both 
P^ and Q^"*"^ are diagonal matrices acting on the space of vectors indexed by {i~ , i+) or 
(*^+i, «^!+i)- On the other hand, the matrices are block diagonal, acting separately 
on the — and + indices. Note that the graphical notation used in (47 1 is not directly 
related to that of sections [2] and ^ . The similarity between the them reflects the 
fact that one can introduce a graphical notation to represent any product of tensors 
(functions of multiple arguments) contracted over some indices (summed or integrated 
over some, possibly repeated, arguments). 

To implement the above algorithm, it is important to compute the precise range 
of the summations, the size of each Mf. matrix, that is, the allowed ranges of 
the if spins, and the ranges of the rip ^ summations in the definitions of P and Q 
in (39 1. Whenever the arguments of either the theta or tetrahedral spin networks 



fail to satisfy certain conditions, these networks evaluate to zero. Therefore, the m , 
if, and Up^q ranges are taken to be the largest such that all necessary conditions are 
satisfied. These conditions are 



and 



0{a, b, c) : tri(a, b, c), 
a b e 



(48) 



tri(c, d, /), tri(a,6, /), tri(a,d, e), tri(c, 6, e), (49) 



[c d f\ 

where the abbreviation stands for the triangle inequality and parity constraint 

tri(a, b,c) : a<b + c,b<c+a,c<a + b, and a + b + c = (mod 2). (50) 
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It can be shown that these conditions, once collected from equations (39 1, (45 1 and 



(46 1, are sufficient to make all summations involved in the algorithm finite. The 
linearity of the triangle inequalities also implies that the upper bound on all sums 
grows linearly with the magnitude of the input i-, j-, and fc-spins. However, it is 
important for an efficient implementation to obtain the tightest possible bounds on 
each of the summation indices. 

The dimension of each M^, is 0{j^) x O(j^), implying d — 2. However, each Me 
decomposes into sparse (diagonal or block diagonal) factors. The filling complexity 
parameters for the largest of these factors, P and Q, are f — 1, F = 2, and F+ f — 3. 
Also, the cost of a matrix-vector multiply is parametrized by Z? = 2 + 1 = 3, giving 
2 + d + D = 7. Therefore, the run time complexity of evaluating an EPR or FK vertex 
amplitude is 0{f). This estimate compares favorably to simply treating Me as a 
dense 0{j^) x 0{p) matrix, which would imply an overall 0{j^) run time complexity. 

5.3 Algorithms for boundary states 

Contracting a boundary state, as described in section |2.3| with the vertex ampli- 



tude ( 40 ) gives the following partition function 

Z^,^ A^{{je,e,le.kle})-^{{jc,e,ie.kle})- (51) 

A naive approach to the problem of computing Zij, would wrap an algorithm to com- 
pute (as described in the previous section) in as many outer sums as there are 
spins in {jc,e, *e, el- Namely, for the BC model, this would produce a calcula- 
tion of run time complexity 0{i^'^^^) = 0(j^^), with 10 outer spin sums. The EPR 
model would yield 0(j^^^^) — 0{j^^), with 15 outer spin sums, and the FK model 
Q^j7+35^ = 0(j^^), with 35 outer spin sums. Clearly, with the naive approach, these 
problems become intractable. Fortunately, when dealing with a factored state (as 



defined tentatively in section 2.3 and more precisely in the following sections), these 



summations may be absorbed into a redefinition of the M,, matrices of sections |5.2.1| 



and 5.2.2 producing again a generalized CE algorithm: 

Z^= 0tr[M4M3M2MiMo], (52) 



where (f> is still defined by equation (41 1 and the sign factor is necessarily absorbed 
into the M^. This approach is described in the next two sections. 

It is important to note that the allowed spin summation ranges, which determine 
the dimensions of the Me matrices, may be strongly impacted by the presence of 
a finitely supported boundary state. It is convenient for our purposes to keep the 
assumption that, even in the presence of a boundary state, the summation range 
for each spin is still of order 0(j), an assumption justified for the boundary states 



proposed in sections 3.2 and 3.3 (or finitely supported approximations to them). The 
run time complexity will be analyzed only for this case. However, the same analysis 
can be easily performed in other cases, where some of the spin summation ranges are 
significantly different. 
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5.3.1 BC vertex with boundary states 

For the BC model, consider a factored boundary state of the form 

'fi{jc,e})=l[i^cAjc,e). (53) 



The dependence of the matrices given in (42 1 on {jc.e} allows us to obtain the form ( 52 1 
with the following redefinition: 



(^e + 1) V^2,e-l(j2,e-l) 
^'(j2,e, ie+l,rn-) 9{j2,e, Je+l, "1+) 



ie J2,e rn 




ie i2,e m+ 


ie+1 J2,e-1 jl,e 




ie+1 j2,e-l jl,e 



31,1, 



0{i2,e-l,ie+l,jl.e) 0{j2,e,ie, jl,e) 

Graphically we represent the above equation as 

ie+l 

Me = 

32,e 



(54) 





Me°"« 











(55) 



where M°"^ corresponds to the right hand side of equation (42 1 and tp refer to the 



appropriate factors of the boundary state ( 53 1 . The tadpole ^ shows an internal 
summation over ji e necessary to form the matrix elements of Me. It is shown here 
to highlight the location of the extra summation insertion and the possible relation of 
"01, e to other spins. Note that, without any modification to the evaluation algorithm, 
we can generalize the notion of factored boundary states to include factors of the form 

■0(jl,e, J2,e, j2,e-l)- 

Notice that in this case is dense and of size 0{p) x 0{j^). Hence, the algo- 
rithm's runtime complexity is 0{j^), as c? = 2, _D = 2 + 2, and 2 + d + D = 8, while 
the filling parameters are f — I and F — A. Interestingly enough, the tets satisfy 
an identity which allows us to decompose Afe into sparse factors speeding up both 
the product-trace and matrix filling, thus reducing the run time complexity to O(j^). 



This identity is known as the Biedenharn- Elliot identity 13 24 
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The product of two tets in equation (|54|) can be rewritten using this identity as 



(56) 
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Hence, we can factor as follows: 



V2,e».+l _ (^e + 1) V'2,e-l02,e-l) V;^ ^/-i^e (jl,e) A, 
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Graphically, this rewriting can be show to be a factorization: 

ie+l - 
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where the factors are given explicitly by 
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The decomposition is not completely unique; some of the terms may be distributed 
differently among the factors. However, in this factorization, the dependence of ipi^e 
can only be generalized to (ji,e, j2,e)- 

Thus Me is clearly decomposed into sparse factors, as each of Ae, Be and Ce is 
dense in some indices, but diagonal in others. Computing the run time complexity, 
we get 0(/), as d = 2, = 2 + 1, and 2 + d + D = 7, while / = 0, F = 3 and 
F + / = 3 for filling either ^ or C. Note that the matrices Bf. contracted with the 
■01, e factors do not depend on m+. Hence, their computation can be done outside the 
m+ summation loops and becomes completely subdominant. 

Curiously, the most practically efficient implementation of the algorithm described 
in this section, as carried out by Christensen 15 , turns out to be a hybrid of 0{j^) and 



0{j'^) versions. The factorization (59 1 greatly speeds up the matrix filling step, while 
the simplicity of the dense version of the product-trace operation is still advantageous 
for all inputs tried to date (up to about jo — 10). 
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5.4 New vertices with boundary states 

For the EPR and FK models, consider respectively 



(63) 



and 



(64) 



Again, we shall only discuss the FK model explicitly, as the EPR model can be directly 
obtained by dropping fc-dependent ips and substituting k = j everywhere else. 



Essentially, we want to compute the quantity Zxj, from equation ( 51 1 with a suitably 



factorable boundary state ^ and the vertex amplitude specified by equation ( 44 ) . This 



expression for Z^j, can be cast into the form ( 52 ) with the following redefinition of 
given directly in graphical form: 



J2, 



"e+l 



Jl.e 



]2.e-l 



Tl 



(65) 



where T denotes a product of T and N from equations ( 45 1 and ( 46 1 . Writing out 



this factorization of with all indices shown explicitly, while straight forward, is 
cumbersome and not particularly enlightening. It should now be clear that, for this 
factorization of the Me, individual factors of the boundary state may depend on clus- 
ters of spins of the form (?e, jl,e, J2,e, fcf e, fcf g) aS well aS (le+l, jl,e, J2,e-1, fc?,e' ^2,e-l); 

which are compatible with possible factorizations of the boundary states proposed in 
sections 13.21 and 13.31 

Each Mf. is of size O(j^) x 0{j'^), hence d = 4. However, because of the sparseness 
of the T", P, and Q factors, each matrix-vector multiply takes O(j^) operations, since 
D = 4 -|- 2 for P and Q multiplies and, equivalently in terms of complexity, £> = 5 -|- 1 
for each T multiply. These numbers are identical for both EPR and FK models. On 
the other hand, filling the P and Q matrices for the EPR model does not involve 
summations over fc-spins. Thus, the EPR filling complexity is parametrized by / = 1, 
P = 5, and F + f = 6, while the FK filling complexity is parametrized by / = 3, 
P = 5, and F + f = 8. The overall runtime complexity of the algorithm is 0{j^^), 
2 + d + D ~ 12, both for the EPR and FK models. By conventional standards, this 
algorithm has a very high polynomial complexity exponent. However, it is still a 
substantial improvement over the naive 0{p^) or 0{j^^) estimates found earlier. 



6 Applications of the algorithms 

The algorithms described in the preceding sections have already been implemented and 
applied in several contexts, other than the results presented in this section. Alesci, 
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Bianchi, Magliaro and Perini [T] have used one variation to extend the original wave 
packet propagation calculations of 29 , both to larger input spins and to different 
kinds of observables (although still keeping the j'-spins frozen). Also, a highly opti- 
mized version of the algorithm presented in section |5.3.1| has been implemented by 



Christensen 15 and used to extract next-to-leading-order asymptotics information 



from the BC graviton propagator (cf. section 3.3 1, as a follow-up to 17. While 



the method used in 17 is capable of handling higher input spins, the advantage of 
the new algorithm is much greater precision, which is better suited for subdominant 
asymptotics analysis. 

Here, we apply the new algorithms to the problems of comparison of amplitude 
asymptotics and of wave packet propagation, described respectively in sections |3.1| 



and 3.2 We have already established that the boundary states proposed in the latter 
section are factored states compatible with the new algorithms. However, being gaus- 
sian, they do not have finite support. Fortunately, strong gaussian decay allows us to 
impose a finite cutoff while maintaining acceptable precision. The cutoff chosen for 
all computations presented below was 2.8 standard deviations about the mean. As a 
consequence, the range of each spin sum involved in the computation is still of order 
0(j), as assumed by our run time complexity estimates. 

6.1 Amplitude asymptotics 



As described in section |3.1| because of the different spin argument structure for each 
of the models under consideration, the comparison of their amplitudes has to be done 



at the level of effective vertex amplitudes given in equation (19 1. The results of the 
maximization procedure outlined previously are given below. 

For the EPR model, we found that the maximum allowed value i — 2j maximizes 
the amplitude. This behavior is illustrated in figure [2] for j = 30. 

For the FK model, we found that, for fixed j and k, the amplitude is again maxi- 
mized by the extreme value i = 2k. See figure |3]for the case j = 30 and k = 15. While 
keeping i at the dominant value 2k, for fixed j, the amplitude is maximized by fc = 1, 
although fc = dominates slightly for very small values of j. This k dependence is 
illustrated in figure |4] for the case j — 30. 

Spin foam quantization is similar in spirit to the discretized path integral approach 
to gravity. As such, the spin foam vertex amplitude is often compared to the gravita- 
tional path integral amplitude: 

A,{j) ^ expiiSn], (66) 

where Sr is the Regge action for gravity evaluated on a discrete geometry described 
by the spins j in the large spin limit. For the BC vertex, this view has turned out to 
be overly simplistic. The relation predicted by careful asymptotic analysis is 

AA.l) - D{j) + ^Ji{3)[exp{^SR) + exp(-i5«)] + • • • , (67) 

where D{j) and /i(j) are non-oscillating functions decaying as and respec- 
tively, with (■■•) representing higher order terms. The dominant asymptotic D(j), 
understood to be due to the contribution of degenerate geometries, masks the desired 
Regge action amplitude [8 11 23 



A natural question is whether the new vertices share the same asymptotic behavior. 
Numerical evaluation of the BC vertex is only sensitive to the dominant asymptotic 
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Figure 2: Effective EPR vertex amplitude: all j — 30, all i equal, satisfying < i < 2j. 
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Figure 3: Effective FK vertex amplitude: all j = 30, all k = f 5, all i equal, satisfying 
< I < 2fc. 
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Figure 4: Effective FK vertex amplitude: all j = 30, all i ~ 2k, all k equal, satisfying 
< A; < j. 
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Figure 5: Large j behavior of the effective vertex amplitudes for the BC, EPR and 
FK models. 
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contribution D{j). The subdominant oscillating Regge action term would become 
important only if D(j) is subtracted or if the vertex amplitude is averaged against 
another oscillatory function, in phase with one of the Regge action terms, as done 



in the graviton propagator calculations 17 34 . While analytical asymptotics for 



the EPR and FK vertices are still missing, we can straightforwardly compare the 
numerical asymptotics of the dominant effective vertex amplitudes of all three models. 
This comparison is made in figure [s] For all models, the data shows no oscillations, 
which means that we are most likely seeing only the D(j) asymptotic term. Note 
that the power laws shown in the figure will change if the edge or face amplitudes 
given in section [2] are modified by j-dependent factors. Such modifications can come, 
for instance, from different choices of exponent k for factors of the form {jf + 1)'^, 
contributing to face amplitudes, as was considered in equation (3) of [l^. Different 
choices of this k would correspond to different choices in the path integral analog of 
spin foams. 

6.2 Wave packet propagation 

As an immediate improvement over previous work, the algorithms presented in this 



paper allow us to show the effect of introducing a non-zero r in ( 23 1 and to compare 



with the calculations of 29], which kept r = 0, freezing all j-spins at the background 



value jo- According to equation (25), the size of r is inversely proportional to the 
parameter a. Figure [6] compares the reference wave packet [cf. (24l] with several 



propagated wave packets (j) (each with a different value of a) depending on the single 
fixed i-spin. The wave packets have been normalized such that their absolute values 
squared sum to 1. The wave packet with the largest value of a is essentially identical 
to the one obtained with all j-spins frozen at jo- In that shown previously 

in [29], the reference state ip resembles the propagated wave packet in shape and mean. 
Unfortunately, as the width of the gaussian factors associated to j-spins increases (a 
decreases), the propagated wave packet quickly departs from -0 in both shape and 
mean. Notably, the mean shifts to a significantly higher value of i. 

Second, we can compare the wave packets propagated by the three different models 
in the 9-1 geometry. Figure [t] shows the reference wave packets ■0 [cf. ( pS] )] and 
propagated wave packets 0, depending on the fixed j-spin and for two choices of a. 
These wave packets are also normalized. While he propagated wave packets seem to 
retain a roughly peaked shape, their mean and width mostly differ significantly from 
the reference wave packet. The only exception is the BC model at a = 5. However, 
this is likely a coincidence that does not appear for other choices of parameter values. 

Lastly, we compare the wave packets propagated by the three different models in 
the 4-6 geometry. In general, the propagated wave packet will depend on the four 
fixed j-spins. Unfortunately, it is impractical to either compute or display functions 
on a 4-dimensional domain. Thus, all calculations have been done with the four fixed 
j-spins set equal. Figure [sjshows the reference wave packet [cf. (23l] and propagated 



wave packets (f>, depending on the common value of the fixed j-spins for two choices 
of a. These wave packets are again normalized. As is clearly seen from the figure, 
the propagated wave packets have in general very little similarity with the reference 
one. More pathological behavior is observed in the FK and BC models, the latter at 
a = 1/2, since none of these curves resemble a well formed gaussian wave packet. In 
all cases, the propagated wave packet has little in common with the reference one. 
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Figure 6: EPR 4-1 propagated {(j)) and reference {ip) wave packets, with jo = 3. 




Figure 7: 9-1 propagated and reference wave packets for different models, with jo = 3. 
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Figure 8: 4-6 propagated and reference wave packets for different models, with jq = 4. 



7 Conclusion and outlook 

We have presented three spin foam models in a unified framework: the standard 
Barrett-Crane (BC) model and two more recent proposals, the Engle-Pereira-Rovelli 
(EPR) and Freidel-Krasnov (FK) models. Their vertex amplitudes were simplified 
and explicitly evaluated using spin network recoupling techniques. Despite the differ- 
ent spin argument structure, we have proposed uniform methodologies for comparing 
the results of calculations of asymptotics of effective vertex amplitudes and of expec- 
tation values of specific spin foam observables among the different models. Building 
on the past work of Christensen & Egan, a fast numerical evaluation algorithm for 
the new vertex amplitudes, both in the absence and in the presence of (factored) 
boundary states, was developed and implemented. The run time complexity of these 
algorithms has been analyzed and shown to be orders of magnitude superior to more 
naive approaches. These algorithms were applied to the problems of extracting the 
asymptotic behavior of the new vertex amplitudes and to checking the behavior of 
propagated semiclassical wave packets. 

The results presented in section |6.1| show that the dominant asymptotic behavior 
of the new vertex amplitudes is non-oscillatory and displays power-law decay very 
similar to the BC model, suggesting dominance of degenerate spin foams as for the 
BC model itself. The power law exponents are estimated in figure |5] It should be 
interesting to explore the asymptotics of these amplitudes with analytical techniques 



as well. It is likely that they will reveal structure similar to (67) 



Although the boundary state proposed in section 2.3 is not ideal (see [35] for a more 
realistic proposal), it has the advantage of belonging to the class of factored states, 
allowing efficient numerical computation with algorithms of section [5] Moreover, the 
proposed state can still be used to gauge the qualitative behavior of unfreezing the 



j-spins. Reference 29 had put forward the conjecture that semiclassical wave pack- 



ets, propagated using the EPR dual vertex amplitude, approximate a certain reference 
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gaussian shape, which was demonstrated under somewhat restrictive conditions. Ap- 
phcation of the numerical algorithms described above allowed a broader investigation 
of this question. The class of factored boundary states encompasses the states used 
in the calculations that gave rise to this conjecture. The results of section |6] indicate 
that a generic factored boundary state does not exhibit good semiclassical behavior. 
It is likely that the presence of correlations between spins, which are not not captured 
by a factored state, improves the agreement between propagated and reference wave 
packets. However, until the use of more complicated states is possible, the evidence 
for the conjecture of Magliaro, Perini & Rovelli remains inconclusive. 

These algorithms have also already been implemented and applied by other au- 
thors, as discussed in section |6.2[ While, several wave packet propagation geometries 
have been examined, there are many other ones. Important questions remain: Is any 
one of them theoretically preferable to the others? What is the impact of directly 
using the unfactorable boundary state of 35 ? Another immediate possibility for fur- 
ther investigation is the computation of the graviton propagator matrix elements in 
the EPR and FK models. 
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